## Merging Recipe 

by.age <- ddply(merged2, # data frame to use
                 "demage", # variable
                 summarise, # function to use
                 mean = mean(police.trust))

anov1 <- lm(mean.rrs ~ white.black + poor.rich + white.black*poor.rich, data = white.black.only)
anova(anov1)

t.test(mv04$mean.govtrust, mv14$mean.govtrust)

str(mv04$mean.govtrust)
str(mv14$mean.govtrust)
cohen.d(mv04$mean.govtrust, mv14$mean.govtrust, na.rm = TRUE)
t.test(blacks.04$dwpimp, blacks.14$dwpimp)
cohen.d(blacks.04$dwpimp, blacks.14$dwpimp, na.rm = TRUE)

t.test

str(blacks.14)

plot(by.age, type = "l", lwd = 2, col = "navyblue", bty = "l")

## 2004 Prep 

setwd("C:/Users/Nathen Huang/Dropbox/QMSS Thesis/Martha's Vineyard/Data/Ready for Analysis")
mv04 <- read.csv("2004 Martha's Vineyard - Cleaned.csv")
mv04$year <- rep(2004)

## Highly specific categorization of individuals as white 
mv04$white.check <- ifelse(mv04$white == 1 & (mv04$blackam != 1 & mv04$blaknoam != 1 & mv04$latino != 1 & mv04$asian != 1 & mv04$native  != 1 & mv04$pacisl != 1	& mv04$other != 1), 1, 0)

## Highly specific categorization of individuals as black
mv04$black.check <- ifelse((mv04$blackam == 1 | mv04$blaknoam == 1) & (mv04$white != 1 & mv04$latino != 1 & mv04$asian != 1 & mv04$native  != 1 & mv04$pacisl != 1 & mv04$other != 1), 1, 0) 

mv04$dincome <- as.numeric(mv04$dincome)
mv04$dempolt <- as.factor(mv04$dempolt) ## coercing political orientation into number

mv04$poor <- ifelse(mv04$dincome < 4, 1, 0)
mv04$rich <- ifelse(mv04$dincome > 6, 1, 0)

mv04$type <- ifelse(mv04$white.check == 1, ifelse(mv04$poor == 1, "poor.white", ifelse(mv04$rich == 1, "rich.white", NA)), ## TRUE STATEMENT
                    ifelse(mv04$black.check == 1, ifelse(mv04$poor == 1, "poor.black", ifelse(mv04$rich == 1, "rich.black", NA)), NA))

# Government/Police

mv04$lgencor1 <- as.numeric(mv04$lgencor1)
mv04$lsupcor2 <- as.numeric(mv04$lsupcor2)
mv04$lpolice3 <- as.numeric(mv04$lpolice3)
mv04$llocpol4 <- as.numeric(mv04$llocpol4)
mv04$lcong7 <- as.numeric(mv04$lcong7)
mv04$llocgov8 <- as.numeric(mv04$llocgov8)
mv04$lpres9 <- as.numeric(mv04$lpres9)

police <- c("lpolice3", "llocpol4")
mv04$police.trust <- rowMeans(mv04[,police], na.rm = TRUE)

legal <- c("lgencor1", "lsupcor2", "lpolice3", "llocpol4", "lcong7", "llocgov8", "lpres9")

gov <- c("lgencor1", "lsupcor2", "lcong7", "llocgov8", "lpres9")
mv04$mean.govtrust <- rowMeans(mv04[,gov], na.rm = TRUE)


# Racial Rejection Sensitivity **** 

rrs <- c("parestcn", "parestex", "papolcn", "papolex", "pashopcn", "pashopex", "paatmcn" , "paatmex")

mv04$rrs1 <- as.numeric(mv04$parestcn * mv04$parestex) 
mv04$rrs2 <- as.numeric(mv04$papolcn * mv04$papolex)
mv04$rrs3 <- as.numeric(mv04$pashopcn * mv04$pashopex)
mv04$rrs4 <- as.numeric(mv04$paatmcn * mv04$paatmex)

rrs4 <- c("rrs1", "rrs2", "rrs3", "rrs4")
mv04$mean.rrs <- rowMeans(mv04[,rrs4], na.rm = TRUE)

## 2014 Prep 

setwd("C:/Users/Nathen Huang/Dropbox/QMSS Thesis/Martha's Vineyard/Data/Ready for Analysis")
mv14 <- read.csv("MV 2014 Data for Analysis.csv")
mv14$X <- NULL
mv14$X.1 <- NULL

mv14$year <- rep(2014)

table(mv14$dedu.10) ## Is the variable lableed as dedu.10? 
table(mv14$year) ## Did I insert the year?

mv14$white.check <- ifelse(mv14$demracecheck == 1, 1, 0)
mv14$black.check <- ifelse(mv14$demracecheck == 2 | mv14$demracecheck == 3, 1, 0)

mv14$dincome <- as.numeric(mv14$dincome)
mv14$dempolt <- as.factor(mv14$dempolt) ## coercing political orientation into number

mv14$poor <- ifelse(mv14$dincome < 4, 1, 0)
mv14$rich <- ifelse(mv14$dincome > 6, 1, 0)

mv14$type <- ifelse(mv14$white.check == 1, ifelse(mv14$poor == 1, "poor.white", ifelse(mv14$rich == 1, "rich.white", NA)), ## TRUE STATEMENT
                    ifelse(mv14$black.check == 1, ifelse(mv14$poor == 1, "poor.black", ifelse(mv14$rich == 1, "rich.black", NA)), NA))


mv14$edcourtesy <- 7 - as.numeric(mv14$edcourtesy)
mv14$edinsult <- 7 - as.numeric(mv14$edinsult)
mv14$edservc <- 7 - as.numeric(mv14$edservc)
mv14$edbettr <- 7 - as.numeric(mv14$edbettr)
mv14$ednsmart <- 7 - as.numeric(mv14$ednsmart)

every.disc <- c("edcourtesy", "edinsult", "edservc", "edbettr", "ednsmart")
mv14$mean.disc <- rowMeans(mv14[,every.disc], na.rm = TRUE)

## Merge 

practice <- mv14$oppolexp
table(practice)

# mv14 range - 1: 35, 2: 50, 3: 74, 4: 101, 5:93 

table(mv04$ophouse)
mv04$ophouse <- 5 - as.numeric(mv04$ophouse)
mv04$oppolice <- 6 - as.numeric(mv04$oppolice)
mv04$oppolexp <- as.numeric(mv04$oppolexp)
mv04$oppolexp <- ifelse(mv04$oppolexp < 6, mv04$oppolexp, NA)
mv04$oppolexp <- 6 - as.numeric(mv04$oppolexp)

opp.disc <- c("ophouse", "oppolice", "oppolexp")
opp.disc <- c("ophouse", "oppolice")
mv04$op.disc <- rowMeans(mv04[,opp.disc], na.rm = TRUE)
mv14$op.disc <- rowMeans(mv14[,opp.disc], na.rm = TRUE)

t.test(blacks.04$op.disc, blacks.14$op.disc)
cohen.d(blacks.04$op.disc, blacks.14$op.disc, na.rm = TRUE)

table(mv14$ophouse)
mv14$ophouse <- 5 - as.numeric(mv14$ophouse)
mv14$oppolice <- 6 - as.numeric(mv14$oppolice)
mv14$oppolexp <- as.numeric(mv14$oppolexp)
mv14$oppolexp <- ifelse(mv14$oppolexp < 6, mv14$oppolexp, NA)
mv14$oppolexp <- 6 - as.numeric(mv14$oppolexp)
mv14$op.disc <- rowMeans(mv14[,opp.disc], na.rm = TRUE)

summary(mv14$op.disc)

mv04$mean.rrs <- as.numeric(mv04$mean.rrs)
str(mv04$mean.rrs)

## Correlating Proxy Measures 
cor(mv04$mean.rrs, mv04$op.disc, use = "na.or.complete")
cor(mv14$mean.disc, mv14$op.disc, use = "na.or.complete")
fix(mv04)


mv04$demage <- as.numeric(mv04$demage)
mv04$dedu.10 <- as.numeric(mv04$dedu.10)
mv04$dmomedu <- as.numeric(mv04$dmomedu)
mv04$ddadedu <- as.numeric(mv04$ddadedu)
mv04$demgend <- as.factor(mv04$demgend)
mv04$demmarg <- as.factor(mv04$demmarg)
mv04$dempolt <- as.factor(mv04$dempolt)
mv04$demskin <- as.numeric(mv04$demskin)
mv04$dnumhse <- as.numeric(mv04$dnumhse)
mv04$dincome <- as.numeric(mv04$dincome)
mv04$demclass <- 6 - as.numeric(mv04$demclass)

mv14$demage <- as.numeric(mv14$demage)
mv14$dedu.10 <- as.numeric(mv14$dedu.10)
mv14$dmomedu <- as.numeric(mv14$dmomedu)
mv14$ddadedu <- as.numeric(mv14$ddadedu)
mv14$demgend <- as.factor(mv14$demgend)
mv14$demmarg <- as.factor(mv14$demmarg)
mv14$dempolt <- as.factor(mv14$dempolt)
mv14$demskin <- as.numeric(mv14$demskin)
mv14$dnumhse <- as.numeric(mv14$dnumhse)
mv14$dincome <- as.numeric(mv14$dincome)
mv14$demclass <- 6 - as.numeric(mv14$demclass)

## Regression Analysis - MOVE THIS THIS IS NOT RELEVANT TO MERGE

lm.pol.04 <- lm(police.trust ~ mean.rrs + demage  + dedu.10 + dmomedu + ddadedu + demgend + dempolt + demskin + dincome + demclass + white.check + black.check, data = mv04)
summary(lm.pol.04)
bptest(lm.pol.04)
coeftest(lm.pol.04, vcov=sandwich)

lm.gov.04 <- lm(mean.govtrust ~ mean.rrs + demage + dedu.10 + dmomedu + ddadedu + demgend + dempolt + demskin + dincome + demclass + white.check + black.check, data = mv04)
summary(lm.gov.04)
bptest(lm.gov.04)

lm.pol.14 <- lm(police.trust ~ perc.disc + demage + dedu.10 + dmomedu + ddadedu + demgend + dempolt + demskin + dincome + demclass + white.check + white.check*dincome + black.check, data = mv14)
summary(lm.pol.14)
bptest(lm.pol.14)
coeftest(lm.pol.04, vcov=sandwich)

library(sandwich)    
coeftest(lm.pol.14, vcov=sandwich)

lm.gov.14 <- lm(mean.govtrust ~ perc.disc + demage + dedu.10 + dmomedu + ddadedu + demgend + dempolt + demskin + dincome + demclass + white.check + black.check, data = mv14)
summary(lm.gov.14)
bptest(lm.gov.14)

t.test(whites.14$perc.disc, blacks.14$perc.disc)
cohen.d(whites.14$perc.disc, blacks.14$perc.disc, na.rm = TRUE)

## Actual Merge

merged.0414 <- merge(mv04, mv14, all=TRUE, sort=FALSE)
colnames(merged.0414)
summary(merged.0414$year)
table(merged.0414$year)

str(merged.0414$white.black)
summary(merged.0414$white.check)
summary(merged.0414$black.check)

white.black.only.0414 <- filter(merged.0414, white.check == 1 | black.check == 1)
table(white.black.only.0414$race.check)

white.black.only.0414$white.black <- ifelse(white.black.only.0414$race.check == 1, "White", "Black")
table(white.black.only.0414$white.black)
table(white.black.only.0414$year)

white.black.only.0414$year <- as.factor(white.black.only.0414$year)

lm11 <- lm(lgencor1 ~ white.black*year, data = white.black.only.0414)
anova(lm11)
aov11 <- aov(lm11)
TukeyHSD(aov11)

lm12 <- lm(lsupcor2 ~ white.black*year, data = white.black.only.0414)
anova(lm12)
aov12 <- aov(lm12)
TukeyHSD(aov12)

lm13 <- lm(lpolice3 ~ white.black*year, data = white.black.only.0414)
anova(lm13)
aov13 <- aov(lm13)
TukeyHSD(aov13)

lm14 <- lm(llocpol4 ~ white.black*year, data = white.black.only.0414)
anova(lm14)
aov14 <- aov(lm14)
TukeyHSD(aov14)

lm17 <- lm(lcong7 ~ white.black*year, data = white.black.only.0414)
anova(lm17)
aov17 <- aov(lm17)
TukeyHSD(aov17)

lm18 <- lm(llocgov8 ~ white.black*year, data = white.black.only.0414)
anova(lm18)
aov18 <- aov(lm18)
TukeyHSD(aov18)

lm19 <- lm(lpres9 ~ white.black*year, data = white.black.only.0414)
anova(lm19)
aov19 <- aov(lm19)
TukeyHSD(aov19)

str(white.black.only.0414$year)
lgencorc <- summarySE(white.black.only.0414, measurevar="lgencor1", groupvars=c("white.black","year"), na.rm = TRUE)
lgencorc <- lgencorc[-c(3,6),]
lpresc$row.names <- NULL

ggplot(lgencorc, aes(x=year, y=lgencor1, fill=white.black)) + 
  geom_bar(position=position_dodge(), stat="identity") +
  geom_errorbar(aes(ymin=lgencor1-se, ymax=lgencor1+se),
                width=.2,                    # Width of the error bars
                position=position_dodge(.9)) +
  scale_fill_manual(values = c("black","bisque","black","bisque"),
                    name="Race of Respondent",
                    breaks=c("Black", "White"),
                    labels=c("Black", "White")) +
  xlab("Year") +
  ylab("Favorability") +
  ggtitle("Perceptions of Court System Favorability") +
  scale_y_continuous(breaks=seq(0, 4, 0.5)) +
  theme_bw()

lpresc.14 <- summarySE(white.black.only.0414, measurevar="lpres9", groupvars=c("white.black","year"), na.rm = TRUE)
lgencorc <- lgencorc[-c(3,6),]
lpresc$row.names <- NULL

ggplot(lpresc.14, aes(x=year, y=lpres9, fill=white.black)) + 
  geom_bar(position=position_dodge(), stat="identity") +
  geom_errorbar(aes(ymin=lpres9-se, ymax=lpres9+se),
                width=.2,                    # Width of the error bars
                position=position_dodge(.9)) +
  scale_fill_manual(values = c("black","bisque","black","bisque"),
                    name="Race of Respondent",
                    breaks=c("Black", "White"),
                    labels=c("Black", "White")) +
  xlab("Year") +
  ylab("Favorability") +
  ggtitle("Perceptions of President Favorability") +
  scale_y_continuous(breaks=seq(0, 4, 0.5)) +
  theme_bw()


ggplot(lpolicec, aes(x=poor.rich, y=lpolice3, fill=white.black)) + 
  geom_bar(position=position_dodge(), stat="identity") +
  geom_errorbar(aes(ymin=lpolice3-se, ymax=lpolice3+se),
                width=.2,                    # Width of the error bars
                position=position_dodge(.9)) +
  scale_fill_manual(values = c("black","bisque","black","bisque"),
                    name="Race of Respondent",
                    breaks=c("black", "white"),
                    labels=c("Black", "White")) +
  xlab("Class Categories") +
  ylab("Favorability of National Police") +
  ggtitle("Favorability of National Police, 2014") +
  scale_y_continuous(breaks=seq(0, 3.5, 0.5)) +
  theme_bw()

which(colnames(merged.0414)=="lgencor1")
which(colnames(merged.0414)=="lpres9")

outcomes <- na.approx(merged.0414[,4:12])

outcomes_fa <- factanal(outcomes, factors = 5, rotation = "varimax")
summary(outcomes_fa)

outcomes_fa


colnames(merge[,46:97])


table(merge$year)

merged.0414$lgencor1 <- as.numeric(merged.0414$lgencor1)
merged.0414$lsupcor2 <- as.numeric(merged.0414$lsupcor2)
merged.0414$lpolice3 <- as.numeric(merged.0414$lpolice3)
merged.0414$llocpol4 <- as.numeric(merged.0414$llocpol4)
merged.0414$lcong7 <- as.numeric(merged.0414$lcong7)
merged.0414$llocgov8 <- as.numeric(merged.0414$llocgov8)
merged.0414$lpres9 <- as.numeric(merged.0414$lpres9)
merged.0414$lsecur11 <- as.numeric(merged.0414$lsecur11)


merged.0414$loppolexp <- as.numeric(merged.0414$oppolexp)
merged.0414$dwpimp <- as.numeric(merged.0414$dwpimp)
merged.0414$dwpphys <- as.numeric(merged.0414$dwpphys)
merged.0414$dwpgendr <- as.numeric(merged.0414$dwpgendr)
merged.0414$dwpcmplx <- as.numeric(merged.0414$dwpcmplx)
merged.0414$dwpmultb <- as.numeric(merged.0414$dwpmultb)
merged.0414$dwpmultw <- as.numeric(merged.0414$dwpmultw)
merged.0414$bgfips <- as.numeric(merged.0414$bgfips)
merged.0414$year <- as.factor(merged.0414$year)

## Clean up the number of variables

merged.0414$RA <- NULL 
merged.0414$RA2 <- NULL  
merged.0414$location  <- NULL
merged.0414$Survey <- NULL


## Maybe don't nullify anything to do with race 
merged.0414$demrace <- NULL   
merged.0414$white <-	as.numeric(merged.0414$white)
merged.0414$blackam <-	as.numeric(merged.0414$blackam)
merged.0414$blaknoam	<- as.numeric(merged.0414$blaknoam)
merged.0414$latino	<- NULL
merged.0414$asian <- NULL
merged.0414$native <- NULL	
merged.0414$pacisl	<- NULL
merged.0414$other <-	NULL
merged.0414$other.description <- NULL	
## End of Possible race vars

merged.0414$demskin <- as.numeric(merged.0414$demskin)	
merged.0414$demage <- as.numeric(merged.0414$demage)
merged.0414$demgend <- as.factor(merged.0414$demgend)
merged.0414$demmarg <-	NULL
merged.0414$demcont <-	NULL
merged.0414$demcit.7	<- NULL
merged.0414$dembirth <- NULL	
merged.0414$demclass <-	6 - as.numeric(merged.0414$demclass) ## Lower number, lower class
merged.0414$dempolt <-	as.factor(merged.0414$dempolt)
merged.0414$dlibcnsv <- as.factor(merged.0414$dlibcnsv)	
merged.0414$dedu.10 <-	as.numeric(merged.0414$dedu.10)
merged.0414$dmomedu <-	as.numeric(merged.0414$dmomedu)
merged.0414$ddadedu <-	as.numeric(merged.0414$ddadedu)
merged.0414$doccup <- NULL
merged.0414$dincome <-	as.numeric(merged.0414$dincome)
merged.0414$dstay <-	NULL
merged.0414$dlocate <-	NULL
merged.0414$dyears <- as.numeric(merged.0414$dyears)	
merged.0414$dparvist <- as.numeric(merged.0414$dparvist) 
merged.0414$drentown <- as.factor(merged.0414$drentown)
merged.0414$dgrndpar <- as.factor(merged.0414$dgrndpar)
merged.0414$census <- NULL
merged.0414$red.blue	<- NULL
merged.0414$cenprob <-	NULL
merged.0414$zipcode <- NULL
merged.0414$RARace	<- NULL
merged.0414$RAskin <- as.numeric(merged.0414$RAskin)
merged.0414$outuszip <- NULL
merged.0414$dlnightlodg <- NULL
merged.0414$pskincol <- NULL
merged.0414$zip <- NULL
merged.0414$crossst <- NULL
merged.0414$dlodgppl <- NULL  
merged.0414$rainitials	<- NULL
merged.0414$date.1 <- NULL

## Possible removal?

merged <- merge

merged$doccup <- as.character(merged$doccup)
str(merged$doccup)

keep <- c("subject.id", "doccup", "mean.rrs", "mean.disc")

merging <- merged[,keep] ## Variables worth keeping only

merged <- merged[,1:46] ## First 46 variables of merge

merged <- merge(merged,merging) ## Attaching nontransfer vars onto 46 vars
table(merged2$year)

summary(mv14$police.trust)
lm.pol.time <- lm(police.trust ~ as.factor(year) + op.disc + demage  + dedu.10 + dmomedu + ddadedu + demgend + dempolt + demskin + dincome + demclass + white.check + black.check, data = merged.0414)
summary(lm.pol.time)
coeftest(lm.pol.time, vcov=sandwich)
bptest(lm.pol.time)

cohen.d(blacks.04$dwpimp, blacks.14$dwpimp, na.rm =  TRUE)

lm.gov.time <- lm(mean.govtrust ~ as.factor(year) + op.disc + demage  + dedu.10 + dmomedu + ddadedu + demgend + dempolt + demskin + dincome + demclass + white.check + black.check, data = merged.0414)
summary(lm.gov.time)
coeftest(lm.gov.time, vcov=sandwich)
bptest(lm.gov.time)

by.year <- ddply(merged2, # data frame to use
                 "year", # variable
                 summarise, # function to use
                 mean = mean(police.trust))

g_by.year <- ggplot(by.year, aes(x = year, y = police.trust)) + geom_line(color = "navyblue")
g_by.year

lmyear <- lm(police.trust ~ year, data = merged2)
summary(lmyear)

library(lmtest)
bptest(lm1)

resettest(lmyear)


# Opportunities in society
merge$opeq <- 5 - as.numeric(merge$opeq) ## Higher value, strongly agree

merge$opdscem <- as.numeric(merge$opdscem)
merge$opdscac <- as.numeric(merge$opdscac)
merge$ophouse <- as.numeric(merge$ophouse)
merge$opshop <- as.numeric(merge$ophouse)
merge$oppolice <- as.numeric(merge$ophouse)

# Diversity Work and Politics 

library(QMSS)
merge$dwpsoc <- ReverseThis(merge$dwpsoc) ## Higher values, more exposure to diversity
merge$dwpwork <- ReverseThis(merge$dwpwork) 
merge$dwpdiv <- ReverseThis(merge$dwpdiv) 

merge$dwppol <- ReverseThis(merge$dwppol)
merge$dwppoldv <- ReverseThis(merge$dwppoldv)
merge$dwprace <- ReverseThis(merge$dwprace)

## Neural Network

mv04_inputs <- model.matrix(~ mean.rrs + demage + dedu.10 + dmomedu + ddadedu + demgend + demmarg + dempolt + demskin + dnumhse + dincome + demclass + white.check + black.check, data = mv04)[,-1]

mv04_target <- mv04$police.trust
str(mv04_target)

mv14_inputs <- model.matrix(~ mean.disc + demage + dedu.10 + dmomedu + ddadedu + demgend + demmarg + dempolt + demskin + dnumhse + dincome + demclass + white.check + black.check, data = mv14)[,-1]

mv14_target <- mv14$police.trust

model_mv <- elman(mv04_inputs, mv04_target, size = c(2,2), learnFuncParams = c(0.1), maxit = 500, inputsTest = mv14_inputs, targetsTest = mv14_target, linOut = TRUE, na.omit)
merge$dwpimprs <- ReverseThis(merge$dwpimprs)

# Main questions 

## 1. lgencor1 
library(effsize)

t.test(blacks.04$lgencor1, blacks.14$lgencor1)
cohen.d(blacks.04$lgencor1, blacks.14$lgencor1, na.rm = TRUE)

t.test(poor.blacks.04$lgencor1, poor.blacks.14$lgencor1)
cohen.d(poor.blacks.04$lgencor1, poor.blacks.14$lgencor1, na.rm = TRUE)

t.test(rich.blacks.04$lgencor1, rich.blacks.14$lgencor1)
cohen.d(rich.blacks.04$lgencor1, rich.blacks.14$lgencor1, na.rm = TRUE)


t.test(poor.04$lgencor1, rich.04$lgencor1)
cohen.d(poor.04$lgencor1, rich.04$lgencor1)

t.test(whites.04$lgencor1, whites.14$lgencor1)
cohen.d(whites.04$lgencor1, whites.14$lgencor1, na.rm = TRUE)

t.test(poor.whites.04$lgencor1, poor.whites.14$lgencor1)
t.test(rich.whites.04$lgencor1, rich.whites.14$lgencor1)

t.test(whites.04$lgencor1, blacks.04$lgencor1)
cohen.d(whites.04$lgencor1, blacks.04$lgencor1)

t.test(poor.blacks.04$lgencor1, rich.blacks.04$lgencor1)
cohen.d(poor.blacks.04$lgencor1, rich.blacks.04$lgencor1)

t.test(rich.blacks.04$lgencor1, rich.blacks.14$lgencor1)
cohen.d(rich.blacks.04$lgencor1, rich.blacks.14$lgencor1, na.rm = TRUE)

t.test(poor.blacks.04$lgencor1, poor.blacks.14$lgencor1)
cohen.d(poor.blacks.04$lgencor1, poor.blacks.14$lgencor1)

results1 = lm(lgencor1 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results1))

lm1 <- lm(lgencor1 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm1)
bptest(lm1)
coeftest(lm1, vocv = hccm(lm1)) ## Robust Standard Errors Correction

## 2. lsupcor2

t.test(whites.04$lsupcor2, whites.14$lsupcor2)
cohen.d(whites.04$lsupcor2, whites.14$lsupcor2, na.rm = TRUE)

t.test(blacks.04$lsupcor2, blacks.14$lsupcor2)
cohen.d(blacks.04$lsupcor2, blacks.14$lsupcor2, na.rm = TRUE)

t.test(poor.blacks.04$lsupcor2, poor.blacks.14$lsupcor2)
cohen.d(poor.blacks.04$lsupcor2, poor.blacks.14$lsupcor2, na.rm = TRUE)

t.test(rich.blacks.04$lsupcor2, rich.blacks.14$lsupcor2)
cohen.d(rich.blacks.04$lsupcor2, rich.blacks.14$lsupcor2, na.rm = TRUE)

t.test(poor.04$lsupcor2, rich.04$lsupcor2)
cohen.d(poor.04$lsupcor2, rich.04$lsupcor2)

t.test(whites.04$lsupcor2, blacks.04$lsupcor2)
cohen.d(whites.04$lsupcor2, blacks.04$lsupcor2, na.rm = TRUE)

t.test(poor.blacks.04$lsupcor2, rich.blacks.04$lsupcor2)
cohen.d(poor.blacks.04$lsupcor2, rich.blacks.04$lsupcor2)

t.test(poor.whites.04$lsupcor2, rich.whites.04$lsupcor2)

results2 = lm(lsupcor2 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results2))

lm2 <- lm(lsupcor2 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm2)
bptest(lm2)
coeftest(lm2, vocv = hccm(lm2)) ## Robust Standard Errors Correction


## 3. lpolice3

t.test(whites.04$lpolice3, whites.14$lpolice3)
cohen.d(whites.04$lpolice3, whites.14$lpolice3, na.rm = TRUE)

t.test(blacks.04$lpolice3, blacks.14$lpolice3)
cohen.d(blacks.04$lpolice3, blacks.14$lpolice3, na.rm = TRUE)

t.test(poor.blacks.04$lpolice3, poor.blacks.14$lpolice3)
cohen.d(poor.blacks.04$lpolice3, poor.blacks.14$lpolice3, na.rm = TRUE)

t.test(rich.blacks.04$lpolice3, rich.blacks.14$lpolice3)
cohen.d(rich.blacks.04$lpolice3, rich.blacks.14$lpolice3, na.rm = TRUE)

t.test(blacks.04$lpolice3, blacks.14$lpolice3)
cohen.d(blacks.04$lpolice3, blacks.14$lpolice3, na.rm = TRUE)

t.test(poor.04$lpolice3, rich.04$lpolice3)
cohen.d(poor.04$lpolice3, rich.04$lpolice3)

t.test(whites.04$lpolice3, blacks.04$lpolice3)
cohen.d(whites.04$lpolice3, blacks.04$lpolice3, na.rm = TRUE)

t.test(poor.blacks.04$lpolice3, rich.blacks.04$lpolice3)
cohen.d(poor.blacks.04$lpolice3, rich.blacks.04$lpolice3, na.rm = TRUE)

t.test(poor.whites.04$lpolice3, rich.whites.04$lpolice3)

results3 = lm(lpolice3 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results3))

lm3 <- lm(lpolice3 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm3)
bptest(lm3)
coeftest(lm3, vocv = hccm(lm3)) ## Robust Standard Errors Correction

# 4. llocpol4 

t.test(whites.04$llocpol4, whites.14$llocpol4)
cohen.d(whites.04$llocpol4, whites.14$llocpol4, na.rm = TRUE)

t.test(blacks.04$llocpol4, blacks.14$llocpol4)
cohen.d(blacks.04$llocpol4, blacks.14$llocpol4, na.rm = TRUE)

t.test(poor.blacks.04$llocpol4, poor.blacks.14$llocpol4)
cohen.d(poor.blacks.04$llocpol4, poor.blacks.14$llocpol4, na.rm = TRUE)

t.test(rich.blacks.04$llocpol4, rich.blacks.14$llocpol4)
cohen.d(rich.blacks.04$llocpol4, rich.blacks.14$llocpol4, na.rm = TRUE)



t.test(blacks.04$llocpol4, blacks.14$llocpol4)
cohen.d(blacks.04$llocpol4, blacks.14$llocpol4, na.rm = TRUE)

t.test(poor.blacks.04$llocpol4, poor.blacks.14$llocpol4)
cohen.d(poor.blacks.04$llocpol4, poor.blacks.14$llocpol4, na.rm = TRUE)

t.test(blacks.04$llocpol4, blacks.14$llocpol4)
cohen.d(blacks.04$llocpol4, blacks.14$llocpol4, na.rm = TRUE)

t.test(poor.04$llocpol4, rich.04$llocpol4)
cohen.d(poor.04$llocpol4, rich.04$llocpol4)

t.test(whites.04$lcong7, blacks.04$llocpol4)
cohen.d(whites.04$llocpol4, blacks.04$llocpol4, na.rm = TRUE)

t.test(poor.blacks.04$llocpol4, rich.blacks.04$llocpol4)
cohen.d(poor.blacks.04$llocpol4, rich.blacks.04$llocpol4, na.rm = TRUE)

t.test(poor.whites.04$llocpol4, rich.whites.04$llocpol4)

results4 = lm(llocpol4 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results4))

lm4 <- lm(llocpol4 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm4)
bptest(lm4)
coeftest(lm4, vocv = hccm(lm4)) ## Robust Standard Errors Correction


# 5. lcong7

t.test(whites.04$lcong7, whites.14$lcong7)
cohen.d(whites.04$lcong7, whites.14$lcong7, na.rm = TRUE)

t.test(blacks.04$lcong7, blacks.14$lcong7)
cohen.d(blacks.04$lcong7, blacks.14$lcong7, na.rm = TRUE)

t.test(poor.blacks.04$lcong7, poor.blacks.14$lcong7)
cohen.d(poor.blacks.04$lcong7, poor.blacks.14$lcong7, na.rm = TRUE)

t.test(rich.blacks.04$lcong7, rich.blacks.14$lcong7)
cohen.d(rich.blacks.04$lcong7, rich.blacks.14$lcong7, na.rm = TRUE)


t.test(whites.04$lcong7, blacks.04$lcong7)
cohen.d(whites.04$lcong7, blacks.04$lcong7, na.rm = TRUE)

t.test(poor.blacks.04$lcong7, rich.blacks.04$lcong7)
cohen.d(poor.blacks.04$lcong7, rich.blacks.04$lcong7, na.rm = TRUE)

t.test(poor.whites.04$lcong7, rich.whites.04$lcong7)

results7 = lm(lcong7 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results7))

lm7 <- lm(lcong7 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm7)
bptest(lm7)
coeftest(lm7, vocv = hccm(lm7)) ## Robust Standard Errors Correction

# 6. llocgov 8

t.test(whites.04$llocgov8, whites.14$llocgov8)
cohen.d(whites.04$llocgov8, whites.14$llocgov8, na.rm = TRUE)

t.test(blacks.04$llocgov8, blacks.14$llocgov8)
cohen.d(blacks.04$llocgov8, blacks.14$llocgov8, na.rm = TRUE)

t.test(poor.blacks.04$llocgov8, poor.blacks.14$llocgov8)
cohen.d(poor.blacks.04$llocgov8, poor.blacks.14$llocgov8, na.rm = TRUE)

t.test(rich.blacks.04$llocgov8, rich.blacks.14$llocgov8)
cohen.d(rich.blacks.04$llocgov8, rich.blacks.14$llocgov8, na.rm = TRUE)


t.test(poor.04$llocgov8, rich.04$llocgov8)
cohen.d(poor.04$llocgov8, rich.04$llocgov8)

t.test(whites.04$llocgov8, blacks.04$llocgov8)
cohen.d(whites.04$llocgov8, blacks.04$llocgov8, na.rm = TRUE)

t.test(poor.blacks.04$llocgov8, rich.blacks.04$llocgov8)
cohen.d(poor.blacks.04$llocgov8, rich.blacks.04$llocgov8, na.rm = TRUE)

t.test(poor.whites.04$llocgov8, rich.whites.04$llocgov8)

results8 = lm(llocgov8 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results8))

lm8 <- lm(llocgov8 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm8)
bptest(lm8)
coeftest(lm8, vocv = hccm(lm8)) ## Robust Standard Errors Correction

nrow(blacks.04)
nrow(blacks.14)

# 7. lpres9

t.test(whites.04$lpres9, whites.14$lpres9)
cohen.d(whites.04$lpres9, whites.14$lpres9, na.rm = TRUE)

t.test(blacks.04$lpres9, blacks.14$lpres9)
cohen.d(blacks.04$lpres9, blacks.14$lpres9, na.rm = TRUE)

t.test(poor.blacks.04$lpres9, poor.blacks.14$lpres9)
cohen.d(poor.blacks.04$lpres9, poor.blacks.14$lpres9, na.rm = TRUE)

t.test(rich.blacks.04$lpres9, rich.blacks.14$lpres9)
cohen.d(rich.blacks.04$lpres9, rich.blacks.14$lpres9, na.rm = TRUE)

t.test(blacks.04$lpres9, blacks.14$lpres9)
cohen.d(blacks.04$lpres9, blacks.14$lpres9, na.rm = TRUE)

t.test(poor.04$lpres9, rich.04$lpres9)
cohen.d(poor.04$lpres9, rich.04$lpres9)

t.test(whites.04$lpres9, blacks.04$lpres9)
cohen.d(whites.04$lpres9, blacks.04$lpres9, na.rm = TRUE)

t.test(poor.blacks.04$lpres9, rich.blacks.04$lpres9)
cohen.d(poor.blacks.04$lpres9, rich.blacks.04$lpres9, na.rm = TRUE)

t.test(poor.blacks.04$lpres9, poor.blacks.14$lpres9)
cohen.d(poor.blacks.04$lpres9, poor.blacks.14$lpres9, na.rm = TRUE)

t.test(poor.blacks.04$lpres9, poor.blacks.14$lpres9)
cohen.d(poor.blacks.04$lpres9, poor.blacks.14$lpres9, na.rm = TRUE)

t.test(rich.blacks.04$lpres9, rich.blacks.14$lpres9)
cohen.d(rich.blacks.04$lpres9, rich.blacks.14$lpres9, na.rm = TRUE)

t.test(poor.whites.04$lpres9, rich.whites.04$lpres9)

results9 = lm(lpres9 ~ white.black + poor.rich + white.black*poor.rich, data= white.black.only)
anova(na.omit(results9))

lm9 <- lm(lpres9 ~ mean.rrs + demskin + race.check*dincome + demage + demgend + dempolt + dedu.10, data = mv04)
summary(lm9)
bptest(lm9)
coeftest(lm9, vocv = hccm(lm9)) ## Robust Standard Errors Correction